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Based on the recently developed picture of an electronic ideal relativistic fluid at the Dirac point, we present 
an analytical model for the conductivity in graphene that is able to describe the linear dependence on the 
carrier density and the existence of a minimum conductivity. The model treats impurities as submerged 
rigid obstacles, forming a disordered medium through which graphene electrons flow, in close analogy with 
classical fluid dynamics. To describe the minimum conductivity, we take into account the additional carrier 
density induced by the impurities in the sample. The model, which predicts the conductivity as a function of 
the impurity fraction of the sample, is supported by extensive simulations for different values of £, the 
dimensionless strength of the electric field, and provides excellent agreement with experimental data. 

Graphene has revealed a series of amazing properties, such as ultra-high electrical conductivity 1 ' 2 , ultra-low 
shear viscosity to entropy ratio 3 , as well as exceptional structural strength, as combined with mechanical 
flexibility 4 and optical transparency 5 . Many of these fascinating properties are due to the fact that, 
consisting of literally one single carbon monolayer, graphene represents the first instance ever of a truly two- 
dimensional material (the "ultimate flatland" 6 ). Moreover, due to the special symmetries of the honeycomb 
lattice, electrons in graphene are shown to behave like an effective Dirac fluid of massless chiral quasi-particles, 
propagating at a Fermi speed of about c ~ 10 6 m/s 3,7 . This configures graphene as an unique, slow- relativistic 
electronic fluid, where many unexpected quantum-electrodynamic phenomena can take place 8 . For instance, 
since electrons are about 300 times slower than photons, their mutual interaction is proportionately enhanced, 
leading to an effective fine- structure constant oc gr = e 2 lhv F ~ 1. As a result of such strong interactions, it has been 
recently proposed that this peculiar 2D graphene electron gas should be characterized by an exceptionally low 
viscosity/entropy ratio (near-perfect fluid), coming close to the famous AdS-CFT lower bound conjectured for 
quantum- chromodynamic fluids, such as quark- gluon plasmas 8 . This spawns the exciting prospect of observing 
electronic pre- turbulence in graphene samples, as first pointed out in Ref.3 and confirmed by recent numerical 
simulations 9 . 

Some of the electrical properties of graphene are still not fully understood, such as the linear increase of the 
electrical conductivity with the number of charge carriers, the existence of a minimum conductivity (see Ref. 10, 
and reference therein), and even the nature of the main scattering mechanism limiting the carrier mobility 11 . In 
fact, classical transport theories, based on short-range scattering of electrons by impurities, predict that the 
electric conductivity in graphene should be independent of the carrier density 12 . Recent works in the field 1314 
have shown that such linear dependence might be potentially explained by treating the impurities as screened 
Coulomb scatterers. Nevertheless, some measurements of the change in the electrical conductivity upon immer- 
sion of graphene samples in high-K; dielectric media differ from this conclusion 1516 . Here, we construct a model for 
describing the electrical conductivity in graphene by using a completely different approach, which is based on the 
recently developed picture of an electronic ideal relativistic fluid at the Dirac point. We demonstrate that, 
although this model is based on a semiclassical theory (it cannot take into account all quantum effects, e.g. 
Landau quantization, quantum hall effects, and quantum interference), it captures the main factors that con- 
tribute to such linear behavior and the appearance of a minimum conductivity. 

Since the most likely relevant limiting factor for the graphene conductivity is still subject of controversy, e.g. it 
can be due to random charged impurity centers 13 or strong neutral defects that induce resonant scattering 1516 , we 
will treat the impurities as hard-spheres, hindering the electron flow (scattering electrons), similarly to the way a 
disordered medium does in the context of fluid dynamics. The choice of hard- spheres is based on the experi- 
mental results by Monteverde et al n 9 which suggested that electrons seem to collide mostly with short range 
scatterers of the size of a few carbon-carbon interatomic distances, like voids, adatoms, etc. Since the relativistic 
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fluid approach is derived from the quantum Boltzmann equation 
(QBE) for graphene 17 , a hydrodynamic description of the conduc- 
tivity can be expected to apply as long as the QBE collision operator 
takes proper account of the Coulomb interactions between electrons. 
Therefore, once Coulomb interactions are included in the viscosity of 
the fluid, the conductivity (which in our case, unlike viscosity, is a 
property of the flow rather than of the fluid) becomes a function of 
the dissipation introduced in the system by the impurities, i.e. the 
electron -impurity interaction. 

Here, we treat graphene as a disordered medium and develop a 
model for its conductivity, as a function of the impurity density 
describing the anomalous dependence of the conductivity on the 
carrier density and the minimum conductivity due to the carrier 
density induced by the impurities. The results are compared with 
experimental data yielding very satisfactory agreement. 

Results 

Electronic fluid in graphene. Our treatment is based on the 
hydrodynamic equations derived by Miiller et al. 3 ' 7 , based on the 
quantum Boltzmann equation for electrons in graphene. This 
analysis delivers the value of the transport coefficients, namely the 
fluid shear viscosity, which is an input parameter in our model. The 
hydrodynamic approach in graphene is valid when the inelastic- 
scattering rate due to electron -electron interactions dominates. 
This is the case at low doping, at high temperatures, and in mode- 
rate fields 18 . However, to neglect electron-phonon interactions, we 
will have to stay at a moderately high temperature of around 100 K 19 . 
In absence of magnetic fields, the quasiparticle distribution function, 
f, evolves according to the quantum Boltzmann equation, 

^+?,-V/, + eg-V s f,= -fl[fl, (1) 

where E is an external electric field, e the electric charge of the 
electron, Q\f s ] a collision operator that takes into account the 

electron -electron interactions, v s = scfc^|fc|, k the wave vector, c 

the Fermi speed (~ 10 6 m/s), and the sign s, not to be confused 
with the entropy density, distinguishes between electrons ( + ) and 
holes ( — ) 18,20 . At equilibrium, the probability distribution function 
becomes the Fermi-Dirac distribution, 

*' *) = ( m \iirr — ' (2) 

V / e (sc\k\-fi)/k B T +l 

where /n is the chemical potential and T denotes the temperature. 
Thus, in the hydrodynamic limit, from Eqs. (1) and (2) one can 
derive the equations for the Dirac electron fluid in graphene: 
dp/dt + V-(pw) = 0, for charge conservation; de/dt + V-[(e+p)u] 
= 0, for energy density conservation and 



Pr 



+ (2-V)2 



_ u dp 



-rjV 2 u = pE, 



(3) 



for momentum conservation. Here, e is the energy density, p the 
pressure, p the charge density, u the drift velocity, p r = (e+p)/c 2 , 
and y\ the shear viscosity. 

For the case of undoped graphene (ft = 0), the presence of charge 
density is due to the thermal energy and can be described by, 



P = Pth = e 



he 



(4) 



However, when there are impurities, they can induce electric poten- 
tials on the graphene sample and a correction due to the chemical 
potential must be considered, 

p = p th d> p (fi/k B T), (5) 

where (f> p is a dimensionless increasing function defined in Ref. 7. 
Note that, in our analytical model, we will use this concept in 



order to introduce a minimum conductivity in the graphene sam- 
ple, where the function <b p will be modeled by a free parameter to 
fit the experimental data and will take into account not only the 
carriers generated by the impurities but also other kind of phe- 
nomena that could contribute to induce carrier density. 
The shear viscosity rj, in Eq. (3), can be calculated using 



M(k B T) 
4Pic 2 a 2 



(6) 



where C n ~ 0(1) is a numerical coefficient, a = e 2 l8hc is the effective 
fine structure constant, 8 the relative dielectric constant of the sub- 
strate, and M the number of species of free massless Dirac particles 3,7 . 
Additionally, the entropy densities can be calculated according to the 
Gibbs-Duhem relation e+p = Ts. These equations have been derived 
under the assumption | u\ < c, and therefore the relativistic correction 
term, oc dp/dt, can be neglected, so that the classical Navier- Stokes 
equations are recovered. Note that, despite the high speed of the 
electrons, |iZ|~0.1c, the Reynolds number remains moderate, due 
to nano-metric size of the samples and the high kinematic viscosity of 
the electronic fluid in graphene. 

Kinematic viscosity. Based on Ref. 3, the dynamic viscosity of 
graphene in a sample of linear size L 0 , is given by Eq. (6). This 
equation can be written in the following form: 



n = C n 



M [k B T\ l h 



n 4a 2 \hCDf J Lq ' 



(7) 



where we have introduced the characteristic frequency Of = c/L 0 , and 
by solving the appropriate quantum Boltzmann equation, it is 
concluded that C n ~ 0.449. Eq. (7) can also be rewritten as 

r\ = C n qj 1 h I / Lq, where qf = hc0fl{k B T). Note that, in order for a 

classical (non quantum) picture of electron fluid to apply, the 
energy of excitations must be much lower than the thermal energy, 
i.e. qf<^l, the so-called collision -dominated regime. Taking a typical 
set of parameters (in MKS units), c= 10 6 , L 0 = 10" 6 , T= 100 K, and 
rj/s~0.2h/k B , we obtain f/~10" 20 . Since the Reynolds number is 
dictated by the kinematic viscosity of the fluid, v, rather than by 
the dynamic one, n=pv, with no need of involving the mass 
density, it is of interest to estimate the kinematic viscosity of the 
electron fluid in graphene. 

To this purpose, we appeal to the definition of the Reynolds num- 
ber as given in Ref. 3, namely: 



Re = 



s/k B k B T u 0 L 0 
rj/h hcof v 0 



(8) 



where v 0 = cL 0 . By writing Re = u 0 L 0 /v and equating with the above, 
we obtain 

hcof s/k B 

v = Vo hfJJ7- (9) 

Using nls = 0.2h/k B and qj ^0.07, we obtain v ~ 10 ~ 2 . To be noted 
that, in spite of its extremely low dynamic viscosity, the kinematic 
viscosity of graphene is about four orders of magnitude larger than 
that of water. These four orders of magnitude are more than com- 
pensated by the large speed of the electrons, which is ultimately 
responsible for the sizeable values of the Reynolds numbers which 
can be achieved in graphene samples at micron scales. For instance, 
by taking u 0 = 0.1 c~ 10 5 m/s, for a sample of 1 micron in length, we 
obtain Re— 20 for the global sample, and about Re— 0.04 on the scale 
of the impurities. 

Analytical model description. In this work, we will treat impurities 
as circular rigid obstacles of diameter d. This choice is not arbitrary, 
but it is based on the fact that some experiments 111516 suggest the the 
main scattering mechanism in graphene is due to strong neutral 
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defects, with a range shorter than the Fermi wavelength, inducing 
resonant scattering. Thus, the diameter d can be interpreted as the 
characteristic length for the range of the interaction electron- 
impurity. 

Let us now assume that the electronic fluid moves in the x dir- 
ection as a consequence of an applied electric field E, and V- (p r u) ~ 0 
(incompressible limit). Therefore, we can calculate the force F acting 

on a single impurity due to the electronic flow, as F = j) Tl-dl, where 

IT is the stress tensor defined by Hy = pd^ + p r UjUj + ny, with the 

viscous tensor = rj(duj/dxj + duj/dxj), and / is a unit vector 
perpendicular to the impurity circle. Here the indices i, j = 1, 2 
denote the coordinates x and y. Thus, solving the equations for the 
Dirac electron fluid, in the steady state (all time derivatives are 
neglected), and following an analogous procedure as in classical 
hydrodynamics 21 , we obtain for the drag force, F D = F x , 

F D = lrjv, (10) 

where X is a dimensionless parameter that depends on the Reynolds 
number. Here, v is the velocity of the fluid very far from the impurity. 
We first simulated single impurities with different sizes and different 
fluid velocities (see Fig. 1), obtaining that a linear approximation is 
justified in the range of relevance to this work. Note that, for a perfect 
fluid (v = 0), there is no drag force. However, from the point of view 
of the kinematic viscosity, electrons in graphene are far from being a 
perfect fluid. Therefore, we expect them to play a crucial effect on the 
drag force controlling the total conductivity of the sample. 

Let us denote by </> the impurity fraction, namely the ratio between 
the area not occupied by the impurities and the total area of the 
system, (j) = 1 — Nnd 2 /(4A), with N the number of impurities in a 
sample of area A. We can thus propose the relation between v and the 
total current density as / = (j>pv. Note that <j>p is the existent electronic 
charge density in the graphene sample, since the volume fraction 1 — 
<j> is occupied by the impurities. 

Let us consider a representative elementary area of the sample of 
length d x in the direction of the flow and d y across it. With J\f circular 
impurities in this area, we can write N ' = 4(1 — (j))d x d y / (nd 2 ). The 
total force acting on the electronic fluid due to impurities (equal to 
the total force acting on the impurities due to the fluid), is given by 
Fm=NFd- Here, the distance between impurities is taken suffi- 
ciently large enough to prevent the flow close to an impurity from 
affecting the flow nearby another impurity. 

In order to describe correctly the physics of graphene, we need to 
include in our model an extra feature. Due to the linear Dirac-Weyl 
spectrum of graphene, and the non-existence of a gap between the 
conduction and valence bands, the slightest amount of impurities or 
an external potential will induce charge carriers in the graphene 
sample 22,23 , see Eq. (5). Thus, in our model, the total amount of 
carriers induced by the impurities will be proportional to the impur- 
ity concentration, (1 — (j))A, leading to an extra carrier density in the 
fraction of the sample occupied by the electronic fluid, </>A.The extra 
carriers are then given by y(l — 0)/y, where y is the proportionality 
constant that characterizes the impurity- electron interaction. 

Summing the forces, namely the Lorentz and drag forces, acting on 
the elementary area leads to: <j)pEd x d y + y(l — (j))Ed x d y / '(j) — F tot = 0, 
and by inserting the value of F tot , we obtain 

pud 1 ( (j) y \ 
F D = -— ~ — t + — r )E. (11) 

Replacing Eq. (10), taking into account that / = (j)pv and Ohm's law, 
we can identify the conductivity as: 

p 2 nd 2 ( 4 2 y\ <j> 2 

a= ^T\—$ + - P ) =°°T^+^> ( 12 ) 

where we have introduced the coefficients <j 0 = p 2 nd 2 l(4r]X), and 
a min = a 0 y/p. This equation represents the key result of our paper. 
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Figure 1 | Drag force acting on a single impurity. Drag force F D acting on 
a single impurity as a function of the graphene flow drift velocity for 
different impurity diameters. The solid lines represent the linear 
dependency of the drag force with the velocity of the fluid. In the inset, the 
dependence of the dimensionless parameter X on the impurity diameter is 
shown. 

An analogous derivation, for fluid dynamics in disordered media, 
can be found in Ref. 24, 25. Note that a 0 also can be written as <7 0 = 
n(e 2 nd 2 c 2 /l2k B TvX), where n = pie is the electronic number density 
and v the kinematic viscosity. From this expression, we can see that 
the conductivity in graphene depends linearly on the carrier density, 
thus implying a constant mobility p = olne, in agreement with 
experimental observations 1 ' 26 . In addition, our model can also 
explain why the mobility remains almost constant in the range of 
temperatures where v oc T" 1 (see Ref. 9), and the presence of a 
minimum conductivity in graphene, second term on the rhs of Eq. 
(12), a min = en^c 2 yl\2k B T\X, which is independent of the carrier 
and impurity densities. Indeed, this model cannot explain, as other 
theoretical models, the sublinear behavior of the conductivity for the 
zero range impurity because, in that range, point defects and bound- 
ary conditions start to be dominant. In addition, it cannot describe 
the electron-phonon interaction either, since these have been 
excluded at the outset. All features above will make the object of 
future extensions of this work. 

Dimensionless numbers. For the numerical validation, and in order 
to obtain general results, we will use dimensionless numbers. For this 
purpose, we can rewrite Eq. (3) alternatively as 

du _ N _> 1 _ u dp n ~ . p - , v 
— + (3-V )U+ -Vp+ —± - -^V 2 u= ^E, 13 
at p r p r c 2 at p r p r 

and therefore we obtain, 

du _ ^ 1 „ u dp p -> . . 

at p r p r c l at p r 

where v represents the kinematic viscosity. Let us define the following 
relations: u = u 0 u\ t = to?, (x, y, z) = L 0 (x\ y\ z'), p r = p r0 p r \ p = 
p 0 p', and E = E 0 E f , where the prime quantities are dimensionless, 
and u 0 , t 0 , L 0 , Pro> Po> and E 0 are characteristic values for the 
respective physical quantities. Thus, using the state equation e = 2p, 

1 2 
we can deduce p = -Pro^p' an( ^ e= ^Pro u o 6 ' '• Replacing these 

relations in Eq. (14), multiplying by — , and using the relation 

Uo 

Uo = — , we obtain, 
to 
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Figure 2 | Speed of the electronic flow. Absolute value of the velocity in 
graphene with multiple impurities, for two different impurity fractions, 
0.952 (bottom) and 0.999 (top). The electric field is applied in the x 
direction (from left to right) and set up to 1.77 V/m. 

±w+\4 *%- ^u>= ^4 (is) 

df 3p r ' r 3 c 2 dt' /u 0 L 0 p r0 fi 2 0 p' r 

To simplify this equation, we can equal the characteristic velocity to 
the Fermi speed, u 0 = c. Therefore, we obtain 

du' _, N 1 „, , W dp' v - PnEoLo p' , . 

— + (u'-V)u' + — -Vp + V 2 W = H0 16 

ot 3p r 3 dt p 0 L 0 p r0 p 2 0 p r 

We can identify two characteristic dimensionless numbers. The first 

HqLq 

one is the Reynolds number, which is, Re = , and the second 

v 

one, which we call "S number" is defined by £ = ^° - 2 ° = 

P r0 U Q Pr0 u 0 

where V 0 = E 0 L 0 is the characteristic electrical potential. Using these 
definitions, we get 

du' 1 . , u' dp' 1 „ n _, „ p' ->. 

dt' v ; 3p r ' F 3 dt' Re p' r 

Note that this equation is dimensionless and therefore the universal 
features of the dynamics of the system are controlled only by the 
numbers Re and 8: the latter measures the strength of the electric 
drive, while the former scales inversely with the dissipation opposing 
this drive. For notational simplicity, we will remove primes, leading 
to 

du 1 „ udp 1 . „ p -> , x 

— + (u-V)u+ — Vp + --f - — V 2 u = £—E, 18 
dt v ; 3p r r 3dt Re p r 

Numerical results. Fig. 2 illustrates the speed of the fluid for two 
different impurity densities, dark and yellow colors denoting low and 
high speeds respectively. An electric field of 1.77 V/m was applied in 
x direction (from left to right). Here we can see that for high impurity 
fraction (see Fig. 2, top), the speed of the fluid presents fluctuations 
on larger scales affecting the total conductivity of the sample. From 
the calculation of the electric current density and the electric field, we 
obtain the Ohm's law, giving a slope a. 

The conductivity a is calculated from the numerical slopes and 
plotted as a function of the impurity fraction. The inset of Fig. 3 
reports the comparison between the analytical solution, using Eqs. 
(12), showing an excellent agreement with the numerical data. 



Figure 3 | Comparison between our analytical model and experimental 
data. Comparison between our results (stars) and experimental data for 
the conductivity cr, as a function of Data from Ref. 1 are represented 
by up and down triangles, from Ref. 23 by circles and squares, and from 
Ref. 3 1 by diamonds and crosses, for electrons and holes respectively. Solid 
lines from bottom to top, theory for separations S = 0 according to Ref. 13, 
our results, and theory for 3 = 0.2 nm according to the previous reference. 
In the inset (top) , we show the conductivity as a function of (f> 2 /( 1 — 4>) , with 
an inset to observe the minimum conductivity. In the inset (bottom) we 
amplify the region close to the Dirac point. 

For the fitting parameters, we obtain er 0 = (9.9 ± 0.1) X 10" 2 e 2 /h, 
and a min = 3.4 ± 0.6. Note that there is a difference between our 
analytical model and the numerical simulations for the minimum 
conductivity. This difference is due to the fact that, for high impurity 
densities, the flow around one impurity starts to affect the flow 
around the others, and therefore, Eq. (12) needs some additional 
terms. In particular, the approximation F tot =Af F D does not hold 
anymore and non-linear correction terms should be included. Thus, 
while the minimum conductivity given by the analytical model is ~ 
3.4 e 2 /h, the simulations give ~ 4 e 2 /h. We have verified that the 
conductivity of graphene, as computed in our model, does not show 
any appreciable dependence on the size of the system. 

In order to compare with experiments, we express the conductivity 
in terms of the ratio nln b where n { = (2.91 X 10 16 m" 2 )(l —(/)), in our 
case. According to this expression and setting n = n 0 , we obtain that 
% = (l/a)«/« f = (l-</>)~\ with a = 4.85 X 10" 3 . Inserting this result 
into Eq. (12), we obtain ala Q = £(1 — 1/y 2 . Note that for values £ > 
1, i.e. n/rii > a = 4.85 X 10" 3 , this equation tends to 

V-Vo£ + Vmin= Vo min . (19) 

a rii 

This corresponds to the linear dependence obtained by different 
theoretical models for graphene 1314 ' 22 ' 23,27 " 30 . 

In Fig. 3, we see the dependence between the conductivity and the 
ratio n/n b and we clearly observe the prediction for the minimum 
conductivity of our model. The experimental data have been taken 
from Refs. 1, 23, 31, and compared with the results of the present 
work, showing good agreement. In Fig. 3, we also compare with the 
model proposed by Hwang et al 13 , where the impurities are located in 
a plane (substrate) parallel to the layer of graphene, with a separation 
3 between the layers. In Fig. 4, we compare our results with Coulomb 
impurity charges in random phase approximation (RPA) 13 . In the 
RPA model, the Boltzmann transport equation is used with impur- 
ities that are located randomly in the graphene sample. Our model 
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shows good agreement in the slope with the RPA model, however, we 
achieve higher values due to the shift made by the minimum con- 
ductivity. 

Discussion 

We have developed an analytical model which accounts for a 
linear behavior of the conductivity with the electron density n, 
as well as with the ratio nln b in the limit «/«,■ > 4.5 X10" 3 . In 
addition, it can also model the minimum conductivity in graphene 
as a consequence of the carrier density induced by the presence of 
impurities. Our model is based on a hydrodynamic description of 
electron flow in graphene, whereby Coulomb interactions are 
included through the viscosity of the electron fluid, and is valid 
in the collision-dominated regime. In this model, the impurities 
are treated as hard-sphere obstacles submerged on the electronic 
fluid, based on the fact that some experiments 111516 suggest that 
strong short-range neutral scatterers are the main scattering mech- 
anism in graphene. Although this idea and the one about the long- 
range Coulomb scatterers are still object of controversy, the fact 
that the present analytical model can account for the conductivity 
of graphene suggests that indeed the short-range scattering models 
might be appropriate for graphene. 

This work is based on the hydrodynamic description of electrons 
in graphene proposed in Ref. 3, 7, which is a model developed for 
undoped graphene that neglects the electron-impurity and electron- 
phonon interactions. Here we have -extended- this approach by 
adding the electron-impurity interactions through a macroscopic 
porous media approach. Since this approach rests on basic conser- 
vation laws, it is supposedly very robust and independent on the 
validity of an underlying quantum Boltzmann equation, so long 
the microscopic interactions justify the build-up of a macroscopic 
viscosity (no superconductivity or other macroscopic quantum 
effects of that sort). Thus, our model is able to reproduce experi- 
mental results to a satisfactory degree of accuracy. 

For the set of parameters investigated in the present work, linear 
Ohm's law appears to apply throughout. However, based on Ref. 9, 
by increasing the size of the impurities (less screening), non-Ohmic 
behavior could occur, due the onset of pre-turbulent phenomena 
within the graphene sample. It would be very interesting to verify 
such possibility by future experiments, as well as the inclusion of the 
electron-phonon interaction to model both, suspended and sup- 
ported samples, at higher temperatures. 

Methods 

For the simulation, we use the hydrokinetic fluid solver proposed by Mendoza 
et a/. 32 " 34 , adapted to two-dimensional flow in graphene 9 . The simulation was 
implemented on a grid of size 256 X 512 cells, representing a rectangular graphene 
sample of size 1.5 X 3 /mi. We set up samples with a fixed number of impurities 
located randomly on the grid, each impurity covering one cell size, varying (f> between 
0.4 and 0.999. The Dirac-quasiparticle fluid in graphene has a kinematic viscosity v = 
8.57 X 10" 3 m 2 /s (see Ref. 9), and by taking the Fermi speed u 0 = 10 6 m/s as a 
characteristic speed, we obtain a Reynolds number Re = 350. Equating Re for both 
systems, in physical and numerical units, the cell size and time step are fixed to Sx = 
5.86 nm and dt = 5.86 fs. For a given temperature, T 0 = 100 K in our case, we can 
calculate the carrier density induced by the thermal energy with Eq. (4), n 0 = 1.41 1 X 
10 14 m" 2 and therefore, using the approximate relation € = 2nok B To X7 , the energy 
density e = 3.90 x 10" 7 J/m 2 and the density p^ = 5.84 X 10" 19 kg/m 2 . Using the 
electron charge, we obtain the charge density, p 0 = 2.26 X 10" 5 C/m 2 . In numerical 
units, these values correspond to n 0 = 4.845 X 10" 3 , e = 2/3, and p r0 = p 0 = 1.0, 
where the charge of the electron is e = 2.064 X 10 2 . Using the characteristic velocity 
u Q , we can calculate the value of the characteristic current density J 0 = p 0 u Q = 22.6 A/ 
m or Jo = PqUq = 1.0 in physical and numerical units, respectively. On the other hand, 
to obtain realistic values of S, we use a characteristic electric field of E 0 = 4.41 V/m, 
which in numerical units corresponds to E 0 = 10" 6 . In this work, £ takes values from 
10" 5 to 2 X 10" 4 . 

To model the extra carrier density induced by the impurities, as described in Eq. 
(5), we introduce an extra density charge Ap localized on each impurity position. 
Therefore, each impurity contributes a quantity p* = Ap dx 2 IA to the total charge of 
the sample, such that p = p th + p* N = p th + Ap n b where n b denotes the impurity 
density. This linear dependence between p and n f is in qualitative agreement with 
experimental data 23 . We made several simulations for different values of Ap, finding 




n/n. 



Figure 4 | Comparison between analytical approaches. Conductivity a as 
a function of nlrii for different types of scattering models 13 . RPA is the 
conductivity calculated by using a random phase approximation with 
Coulomb scatterers. The unscreened Coulomb interaction would yield a 
conductivity smaller than the minimum value in graphene, over the entire 
range of gate voltages. 

that Ap = 60 leads to a minimum conductivity of 4 e 2 /h. In our analytical model, this 
value corresponds toy = 60, in numerical units. The simulations ran up to 5 X 10 5 
time steps. 
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